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Abstract 

We study the extended Thomas- Fermi (ETF) density functional of the su- 
perfluid unitary Fermi gas. This functional includes a gradient term which is 
essential to describe accurately the surface effects of the system, in particu- 
lar with a small number of atoms, where the Thomas-Fermi (local density) 
approximation fails. We find that our ETF functional gives density profiles 
which are in good agreement with recent Monte Carlo results and also with 
a more sophisticated superfluid density functional based on Bogoliubov-de 
Gennes equations. In addition, by using extended hydrodynamics equations 
of superfluids, we calculate the frequencies of collective surface oscillations of 
the unitary Fermi gas, showing that quadrupole and octupole modes strongly 
depend on the number of trapped atoms. 
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I. INTRODUCTION 



In a system of spin 1/2 interacting fermions, the unitary regime is commonly referred 
to as the situation in which the s-wave scattering length a greatly exceeds the average 
interparticle separation, thus n\a\^ ^ 1, where n is the number density [1]. In 2002 it was 
shown experimentally that such systems are (meta) stable [2], and they have been studied 
extensively ever since [1] . On the basis of very accurate Monte Carlo results, we have recently 
introduced a simple but reliable extended Thomas-Fermi (ETF) density functional for the 
unitary Fermi gas at zero temperature [3]. 

In this paper we show new results obtained with our ETF functional. In particular, we 
analyze the density profiles of the balanced zero-temperature unitary Fermi gas confined by 
a harmonic trap, finding that our ETF functional gives density profiles which are in good 
agreement with recent Monte Carlo results [4,5] and with a more sophisticated superfluid 
density functional based on Bogoliubov-de Gennes equations [6]. We then obtain extended 
hydrodynamics equations of superfiuids by considering small, time-dependent perturbations 
of the ETF functional around its ground-state. We show that, under isotropic harmonic 
confinement, the quadrupole and octupole frequencies of a Fermi gas at unitarity depend on 
the number of confined atoms, contrary to the monopole and dipole frequencies, which are 
solely determined by the external potential. 

II. EXTENDED THOMAS-FERMI DENSITY FUNCTIONAL 

At zero temperature the Thomas-Fermi (TF) energy functional [1] of a dilute and ultra- 
cold unpolarized two-component Fermi gas trapped by an external potential U (r) is 



with e{n;aF) bulk energy per particle, n{r) = n^{r) + nj(r) total density {ni{r) = nj^(r)) 
and ap the s-wave scattering length. The total number of fermions is 




(1) 
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Jdh n{r) . (2) 

By minimizing Etf one finds 

li{n{r);aF) + U{r) ^ jl , (3) 

with fj,{n;aF) — ^^"'^q'^"'^^'* bulk chemical potential of a uniform system and p chemical 
potential of the non uniform system. 

For the uniform unitary Fermi gas [7] the s-wave scattering length ap diverges: 

ap — >■ ±C)0 , (4) 

and the only length characterizing the uniform system is the average distance between 
particles n~^^^. In this case: 

^(n; = 4^(3^^/^ = 4^^, (5) 

with ep Fermi energy of the ideal gas and ^ a universal parameter. The bulk chemical 
potential associated to Eq. (5) is 

The TF functional must be extended to take into account other characteristic lengths 
related to the spatial variations of the density, besides the average particle separation. As 
a consequence, the energy per particle must contain additional terms, which scale as the 
square of the inverse of these various lengths. For this reason, as a simple approximation, 
we add to the energy per particle the term 

8m v? 2m n 

Historically, this term was introduced by von Weizsacker [8] to treat surface effects in nuclei. 
Moreover, according to the Kirzhnits expansion of the quantum kinetic operator in powers 
of h [19], A must take the value A = 1/9 [19,20] for an ideal, noninteracting, Fermi gas. 
Here we consider A as a phenomenological parameter accounting for the increase of kinetic 



energy due the spatial variation of the density. We also observe that the TF functional 
has a pathological behavior at the surface. In fact, the TF density profile becomes zero at 
a finite distance ttf from the center of the cloud, while according to quantum mechanics 
the density must goes to zero only at infinity. This pathology of the TF approximation is 
efficiently cured by the inclusion of the gradient term (7). 

Other recent density-functional methods for unitary Fermi gas have been proposed in 
the last few years: the Bogoliubov-de Gennes (BdG) superfluid local-density approximation 
(SLDA) of Bulgac [6] and the Kohn-Sham (KS) density functional approach of Papenbrock 
[9]. We wish to point out that both the energy functionals proposed by Bulgac [6] and 
Papenbrock [9] are functionals of the density through single particle orbitals (the BdG or KS 
orbitals) . Therefore they can be used in actual numerical calculations only when the number 
of fermions is small, since they require a self consistent calculation of single-particle states 
whose number increases rapidly with the number of particles [17,18]. For completeness, we 
remark that another density functional for the unitary Fermi gas has been proposed on the 
basis of the Haldane-Wu statistics [10,11]. 

Our energy functional, that is the extended Thomas-Fermi (ETF) functional of the uni- 
tary Fermi gas, reads 

[d'r n{r) [A— ^^^^ + ^^^(SnYMrf^' + U{r)] . (8) 

We stress that one encounters no severe limitation in the number of particles which may 
be treated through ETF functionals, since in this case the functional depends only on a 
single function of the coordinate, i.e. the particle density. By minimizing the ETF energy 
functional one gets: 



This is a sort of stationary 3D nonhnear Schrodinger (3D NLS) equation. 

The constants ^ and A should be universal i.e. independent on the confining potential 
U{r). The value of the coefficient A is debated. In the papers of Kim and Zubarev [12] 



and Manini and Salasnich [13] the authors set A = 1 over the full BCS-BEC crossover. 
More recently we have suggested A = 1/4 [14-16]. Note that A ~ 0.25 was the prehminary 
prediction at unitarity of effective field theory [21]. We stress, however, that the more recent 
and accurate prediction of effective field theory based on e-expansion is A ~ 0.17 [22]. 

III. FINDING THE UNIVERSAL PARAMETERS OF THE ETF FUNCTIONAL 

To determine ^ and A we look for the values of the two parameters which lead to the 
best fit of the ground-state energies obtained by Monte Carlo data for the unitary Fermi gas 
in a spherical harmonic potential 

U{r) = ]-muj\'^ . (10) 

We use the more recent and reliable Monte Carlo results with N even (complete superflu- 
idity): the fixed-node diffusion Monte Carlo (FNDMC) of von Stecher, Greene and Blume 
[4] . After a systematic analysis [3] we find 

i = 0.455 and A = 0.13 

as the best fitting parameters in the unitary regime. The value ^ = 0.455 coincides with that 
obtained by Perah, Fieri, and Strinati [23] with beyond-mean-field extended BCS theory. 
Fixing ^ = 0.44 (the Monte Carlo result for a uniform Fermi gas of Carlson et al. [24]) we 
find instead A = 0.18 [3]. 

In Fig. 1 we plot the ground-state energy E of the unitary Fermi gas as a function of the 
number N of atoms under harmonic confinement of frequency uj. In the figure we compare 
the FNDMC data with even N [4] with different density functionals: ETF results with best 
fit {i = 0.455 and A = 0.13); ETF results obtained from ^-expansion [22] (^ = 0.475 and 
A = 0.17); and also the TF (LDA) results (^ = 0.44 and A = 0). The figure clearly shows 
that the simple TF predictions are not fully reliable, while the best-fit ETF functional is 
extremely accurate. 
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In Fig. 2 we report the comparison between the density profiles obtained with our best- 
fit ETF functional and the FNDMC results [4,5]. We have also included the SLDA results 
[6]. The figure shows that, near the surface, both ETF and SLDA methods are in very good 
agreement with the FNDMC predictions. Close to the center there are instead appreciable 
differences. Nevertheless, it is important to remind that in the FNDMC approach the 
statistical uncertainty is quite large at small distances from the center of the trap. Moreover, 
also SLDA values for central density depend quite a bit on the size of the BdC basis. 

For the sake of completeness, in Fig. 3 we plot the zoom of the density profile n(r) 
near the surface for N = 20, including also the TF prediction (dot-dashed line). The figure 
clearly shows that the TF data are not reliable near the surface. 

In our determination of ^ and A we have analyzed the unitary gas with an even number 
N of particles. Monte Carlo calculations show a clear odd-even effect (zig-zag effect): the 
ground state energy of N odd particles in the isotropic harmonic trap is 



where the splitting AE^ is always positive. Dam Thanh Son has suggested [25] that, given 
the superfluid cloud of even particles, the extra particle is localized where the energy gap 
is smallest, which is near the edge of the cloud. Dam Thanh Son has also found [25] that, 
for fermions at unitarity, confined by a harmonic potential with frequency uu, the odd-even 
splitting grows as 



where 7 is an unknown dimensionless constant. After a systematic investigation of the 




(11) 



(12) 



FNDMC data we find that 7 = 0.856 gives the best fit [3]. 



IV. EXTENDED SUPERFLUID HYDRODYNAMICS 



Let us now analyze the effect of the gradient term (7) on the dynamics of the super- 
fluid unitary Fermi gas. At zero temperature the low-energy collective dynamics of this 
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fermionic gas can be described by the equations of extended [26,27] irrotational and inviscid 
hydro dynamics : 

^ + V-(nv) = 0, (13) 

m^v + v[- + nin; + ^(r) + ^v'] = . (14) 

ot L zm \jn 1 J 

They are the simplest extension of the equations of superfiuid hydrodynamics of fermions 
[1] where A = 0. 

The extended hydrodynamics equations can be written in terms of a superfiuid nonlinear 
Schrodinger equation (NLSE), which is Galilei- invariant [3]. In fact, by introducing the 
complex wave function 



^(r,t) = ^n(r,t)e*^("-'*), (15) 

which is normalized to the total number N of superfiuid atoms, and taking into account the 
correct phase-velocity relationship 

v(r,t) = ^V^(r,t), (16) 

where 9{r,t) is the phase of the condensate wavefunction of Cooper pairs, the equation 



d 



'^m^^ -i-^^ "^ + 2Ml^h + (1 - 4A)^^J^ , (17) 

is strictly equivalent to the equations of extended hydrodynamics. 

The extended hydrodynamics equations are the Euler-Lagrange equation of the following 
Lagrangian density 

C = -n(e + ^{Ver + U{r) + e{n- i) + \^^^) , (18) 
^ 8m 8m n-^ ' 

which depends on the total number density n(r, i) and the phase ^(r, f) [3]. In the case 
A = it is called the Popov Lagrangian of superfiuid hydrodynamics [28]. 

Using, as previously, Eq. (15) we find that the extended Popov Lagrangian (18) is 
equivalent to the following one: 
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V. SOUND VELOCITY AND COLLECTIVE MODES 



Prom the equations of superfluid hydrodynamics one finds the dispersion relation of 
low-energy collective modes of the uniform {U{r) = 0) unitary Fermi gas in the form 

i - /I- ■ P") 

where Q is the collective frequency, q is the wave number and 

V. - (21) 
V m 

is the Fermi velocity of a noninteracting Fermi gas. 

The equations of extended superfluid hydrodynamics (or the superfluid NLSE) give also 
a correcting term [3], i.e. 

^ /. ^3A hq 2 .... 

which depends on the ratio A/^. 

In the case of isotropic harmonic conflnement, Eq. (10), we study numerically the col- 
lective modes of the unitary Fermi gas by increasing the number N of atoms. By solving 
the superfluid NLSE we find that the frequency Qq of the monopole mode (/ = 0) and the 
frequency fli dipole mode (/ = 1) do not depend on A^: Qq — 2a; and fli — cu, as predicted 
by Castin [29]. 

We find instead that the frequencies 0,2 and Q3 of quadrupole (/ = 2) and octupole 
(/ = 3) modes depend on and on the choice of the gradient coefficient A. The frequencies 
are computed by solving the time-dependent Schrodinger equation (17). We consider the 
initial state 

V'(r,t = 0) = e-«'^,3(r), (23) 

where V'gs(i") is the ground-state wave function, e is a small parameter and Qi is the operator 
which excites the surface mode with angular momentum / (and mi = 0) [17]. For the 
quadrupole mode one has 
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Q2 = -x^-y^ 



(24) 



while for the octupole mode the operator is 



Qs = z{2z'' - - 3y^) . 



(25) 



Both Q2 and Q3 are expressed in cartesian coordinates [17]. In this way, a rehable value for 
the mode frequency can be obtained by following the real-time evolution of the system for 
rather short times, corresponding to few periods of oscillations at most. 

In Fig. 4 we plot the quadrupole frequency as a function of the number N of atoms 
under harmonic confinement. We consider three values of the gradient coefficient A. For 
A = one finds = \f^'jj [1] . As shown by the figure, this is the asymptotic value one finds 
only for a very large N with a non-zero A. Thus, we predict that finite-number effects can be 
detected in the unitary Fermi gas by measuring the quadrupole mode Q2- This finite- number 
effect is stronger for surface modes with a larger angular quantum number /. This is shown 
in Fig. 5, where we plot the octupole mode Q3 as a function of N . 



We have introduced an extended Thomas- Fermi (ETF) functional for the trapped unitary 
Fermi gas. By fitting FNDMC calculations we have determined the universal parameters ^ 
and A of ETF functional. ETF functional can be used to study ground-state density profiles 
in a generic external potential C^(r). We have also introduced a time-dependent version 
of the ETF functional: the extended superfluid hydrodynamics (or superfluid NLSE). The 
superfiuid NLSE can be used to investigate collective modes also for a small number of atoms. 
Our predictions suggest that surface effects are indeed very important for the unitary Fermi 
gas and that the physical observables strongly depend on them. Next-future experiments 
with less than one hundred degenerate Fermi atoms [30] across a Feschbach resonance can 
surely test our theoretical results. 



VI. CONCLUSIONS 
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Fig. 1. Ground-state energy E for the unitary Fermi gas of atoms under harmonic 
confinement of frequency lo. Symbols: FNDMC data with even [4]; sohd hne: ETF results 
with best fit = 0.455 and A = 0.13); dashed line: ETF results obtained from e-expansion 
[22] {i = 0.475 and A = 0.17); dot-dashed lines TF (IDA) results = 0.44 and A = 0). 
Energy in units of hu. 
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Fig. 2. Unitary Fermi gas under harmonic confinement of frequency u. Density profiles 
n{r) for fermions obtained with our ETF (sohd hues), SLDA [6] (dashed hues), and 



FNDMC [5] (filled circles). Lengths in units of an = Jfi/{muj). 
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Fig. 3. Unitary Fermi gas with = 20 atoms under harmonic confinement of frequency 
uo. Density profile n{r) near the surface obtained with our ETF (sohd fine), SLDA [6] 
(dashed fine), TF (dot-dashed fine), and FNDMC [5] (filled circles). Lengths in units of 
an = Jh/{muj). 
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Fig. 4. Quadrupole frequency Q2 of the unitary Fermi gas = 0.455) with atoms 
under harmonic confinement of frequency u. We use three different values of the gradient 
coefficient A. Notice that for A = (TF limit) one finds f22 = V2lj. 
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Fig. 5. Octupole frequency of the unitary Fermi gas = 0.455) with atoms under 
harmonic confinement of frequency u. We use two different values of the gradient coefficient 
A. Notice that for A = (TF limit) one finds fia = \/3lj. 
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